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We present a quenched lattice QCD calculation of the a and (5 parameters of the proton decay 
matrix element. The simulation is carried out using the Wilson quark action at three values of the 
lattice spacing in the range a k, 0.1—0.064 fm to study the scaling violation effect. We find only mild 
scaling violation when the lattice scale is determined by the nucleon mass. We obtain in the contin- 
uum limit, |a(NDR,2GeV)|=0.0090(09)(t^9)GeV^ and |/3(NDR,2GeV)|=0.0096(09)(l®o)GeV^ with 
a and /3 in a relatively opposite sign, where the first error is statistical and the second is due to the 
uncertainty in the determination of the physical scale. 
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Proton decay (or nucleon decay in general) is a char- 
acteristic consequence of grand unified theories (GUTs) 
because of the unification of quarks and leptons into the 
same gauge multiplet. However, no clear evidence of such 
decay process has been observed up to now in spite of con- 
tinual experimental efforts over several decades. Most re- 
cent experimental lower bound of the lifetime is given by 
the Super-Kamiokande experiment: 4.4 x 10'^^ years for 
p — > e"*" -I- 7r° mode and 1.9 x lO'^^ years for p ^ D -\- 
mode at 90% confidence levelQ. Although some naive 
GUT models is already ruled out by this experimental 
bound, we still have several viable GUTs which allow 
the longer proton lifetime at 0{1Q^^^^'^) 2]. Further im- 
provement of the experimental bound could give strong 
constraints on these GUT models. 

One of the main sources of uncertainties in the theoret- 
ical predictions is the evaluation of the hadronic matrix 
elements for the nucleon decays {PS\0\N), where PS 
and N stand for the pseudoscalar meson and the nucleon, 
respectively, and O is the three-quark operator violating 
the baryon number. The conventional procedure of esti- 
mating the hadronic matrix element is to invoke current 
algebra and PCAC and to reduce the three-body matrix 
element into the two-body transition element (OjOliV), 
leaving aside the question as to the validity of PCAC 
with a long extrapolation. Varieties of models have been 
employed to estimate this transition elements, but results 
vary by an order of magnitude; see . 

A promising method to reduce the uncertainty is to 
resort to lattice QCD, which allows direct evaluation of 
non-perturbative effects and has been successfully used in 
giving various weak interaction matrix elements. There 
are already a few calculations to evaluate the two-body 



transition element and even a few attempts to eval- 

uate directly the three body amplitude {PS\0\N) 6, 73- 
Gavela et al.0 argued that the three-body amplitude 
gives proton decay lifetime that differs largely from the 
one derived from a two-body calculation with the use 
of PCAC. The JLQCD calculation 7], however, showed 
that their results are due to a neglect of one of the two 
form factors and that the three-body and two-body cal- 
culations yield the results that agree at a reasonable ac- 
curacy, say 20-30%. 

Lattice QCD calculations, being carried out today, 
however, contain a number of sources that lead to sys- 
tematic errors, such as the quenching approximation, fi- 
nite lattice spacing, finite lattice size, chiral extrapolation 
and so forth. Particularly worrisome are the finite lattice 
spacing effects that could modify the continuum results 
even by a factor 2 if scaling violation is substantial in the 
relevant quantity. In fact, the recent preliminary result of 
the RBC collaboration y gives the matrix elements that 
differs from those by JLQCD ItI by Ri50%, which urges 
us to study the issue of systematic errors. 

In this paper we focus on the issue of the lattice spacing 
effects, by carrying out simulations at three different val- 
ues of bare coupling constant, adopting the lattices that 
are large enough so that finite lattice effects are negligible 
even for baryons, and borrowing the results of a large- 
scale simulation of CP-PACS collaboration for quenched 
hadrons :9] . We consider two-body matrix elements 

(0|eyfe(u'^CPfld^)Ptu'=b(fc = 0)) - aPLUp, (1) 
{Q\e,,k{v^^CPLd')PLu''\p{k^Q)) = PPlUj,, (2) 

expressed by a and (3 parameters, where i, j and k are 



color indices, C is the charge conjugation matrix, Pr/l 
is chiral projection operator and Up denotes the proton 
spinor with the zero spatial momentum. We deal with 
the two-body matrix elements in view of the feasibility 
on current computers, rather than three-body matrix ele- 
ments, which need three-point correlators with finite spa- 
tial momenta injected to disentangle the relevant and ir- 
relevant form factors 0. 

The continuum operators relevant to the a and (3 pa- 
rameters are connected with the lattice operators as 
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and the mixing operator 
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appears due to explicit chiral symmetry breaking of the 
Wilson quark action. The renormalization constants Z , 
Zmix, and Z'^-^^ are evaluated perturbatively at one- loop 
order H The continuum operators are defined in 

naive dimensional regularization (NDR) with the MS 
subtraction scheme. The matrix elements defined on 
the lattice are converted to those in the continuum at 
fi = 1/a and are evolved to /i=2GeV using the two-loop 
renormalization group in the contimmm|lll|. 

To obtain the matrix elements, we consider the ratio 
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where Jp{x, t) is a local sink operator for the proton with 
spin s and Jp_s(0) is the smeared source, 

S,y,z 

,,fe (u'^{x,t)Cj5dHy,t)) u':{z,t) (8) 
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with the smearing function The factor y^'Z^ defined 
by 



\JpAO,0)\p^''Hk^O))^y/Z~pui'''^ 



(9) 



is obtained from the proton correlator with the local 
source and local sink. It is recognized that the precise 
determination of is not easy because of large sta- 

tistical fluctuations of the local-local correlator (see, e.g., 
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FIG. 1: Time dependence of Rnit) (the factor ^J'Z^ being 
removed) for the heaviest (circle) and lightest (triangle) quark 
masses at /9=6.10. 



Fig. 16 of Ref. [9j). On the other hand, the ratio of two- 
point functions in eq.®, calculated using the smeared- 
local proton correlator, is determined well with small sta- 
tistical errors. 

Under this circumstance we calculate y/Z^ from the 
proton correlators generated in high statistics calcula- 
tions of the quenched light hadron spectrum performed 
by the CP-PACS collaboration^. We then carry out new 
simulations with the same parameters to obtain the ra- 
tio of two-point functions including the mixing operator, 
for which we do not necessarily need very high statistics. 
We attain a few percent statistical accuracy for the lat- 
ter, while the overall accuracy is still limited by the error 



of ^/Zp. 

Our simulation generates quenched gauge configura- 
tions at /3=5.90, 6.10 and 6.25, which correspond to lat- 
tice spacings in the range a ~ 0.1—0.064 fm when de- 
termined from the p meson mass mp=0.7684 GeV. The 
spatial lattice size is kept at about 3 fm to avoid finite 
size effects. We take four quark masses corresponding 
to mps/mv ~ 0.75—0.5 for each (3. These parameters 
are the same as those of the CP-PACS spectrum calcu- 
lations , except that we drop the finest lattice and the 
lightest quark mass at each /3 for the computational cost. 
The simulation parameters are presented in TableQ] The 
number of configurations in our simulation is about 1/3 
that of the CP-PACS calculation. We employ for ^E* in 
eq.® the pion quark wave function, which is measured 
for each hopping parameter on 30 gauge configurations 
fixed to the Coulomb gauge except for the t=0 time slice 
where the wall source is placed [l^- 

To estimate \fZp^ we fit the smeared-local proton cor- 
relator to a single exponential Zpexp{—mpt), and then 
fit the local-local proton correlator to Zpexp(— mpt) with 
nip fixed to the value determined from the smeared- 
local correlator, which is borrowed from the CP-PACS 
simulation 9]. Figure ^ shows the ratio of two-point hmc- 
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TABLE I: Simulation parameters and results. The lattice spacing a[fm] in the third column is determined from mp. Results 
are given with three different input quantities for the lattice spacing, i.e. niN, mp, and 
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input 
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x56 0.1020(8) 


300/800 


0.01026(31) 


0.01265(43) 
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0.01312(44) 


0.01621(76) 


6.10 




x70 0.0777(7) 


200/600 


0.01041(35) 


0.01152(41) 
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0.01467(85) 
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48'^ 


x84 0.0642(7) 
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0.00956(35) 
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0.01388(78) 
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FIG. 2: Chiral extrapolation of the a parameter with linear 
(solid) and quadratic (dotted) functions at (3=6.10. 



tions in eq.© with the -^/^ factor removed for the 
heaviest(ii: = 0.15280) and the hghtest(if = 0.15440) 
quark masses at f3 — 6.10. The horizontal lines repre- 
sent the fits together with one standard deviation errors, 
which are smaller than 1%. In Fig.|3we plot the quark 
mass dependence of the a parameter at (3 = 6.10, which 
is well described by a linear function. A similar quark 
mass dependence is observed for the f3 parameter. We 
find that linear plus quadratic extrapolations yield re- 
sults consistent within error bars in the chiral limit at all 
lattice spacings. The a and P parameters in the chiral 
limit obtained by linear extrapolations are summarized 
in Table The errors are at most a few percent. The 
contribution of the mixing operator in eq.|(2Jl is smaller 
than 10%. 

We present in Fig. O a and p in physical units as a 
function of a. We examine three choices, the nucleon 
mass TOjv, the p mass mp, and the pion decay constant 
fjr to determine the physical scale of the lattice. We take 
these physical parameters given by CP-PACS P, because 
a and /3 have dimension 3 and the error of mass scale is 
magnified by a factor of 3, so that high statistics results 
are essential. Figure 13 indicates that scaling violation a 
and (3 is minimized if nucleon mass is used as input. The 
use of mesonic quantities, uip or on the other hand, 
leads to substantial scaling violation. A simple linear 



extrapolation to the continuum limit results in a and /3 
that vary up to 30% depending on the input physical 
scale as found in Table 

We adopt the a and (3, extrapolated to a = 0, using 
the niN input as our central value, since small scaling 
violation would minimize the error associated with the 
continuum extrapolation, and include the uncertainty in 
the physical scale as systematic error. We obtain 

|a(NDR,2GeV)| = 0.0090(09) (+^g) GeV^ (10) 
|/3(NDR,2GeV)| = 0.0096(09) (t^g) GeV^ (11) 

where the first error is statistical and the second one is 
systematic. Since the CP-PACS spectrum calculation^ 
is superior to this work in controlling the systematic er- 
rors using finer lattices and lighter quark masses than 
this simulation, we estimate the ambiguity due to scale 
setting from their results of quenched light hadron mass 
spectrum. They show that the values of mp and /^r in 
quenched QCD deviate from the experiment by -f 7% and 
—2% respectively in the continuum limit, once we set the 
lattice spacing by m^y. The errors of mass scale in a and 
f3 are magnified by a factor of 3 and found to be compa- 
rable with the variation of the results at the continuum 
limit in Tabled This implies that the systematics from 
the physical scale dependence are mostly ascribed to the 
quenching effects. We note that the sign of a and f3 
are relatively opposite, while the overall sign is a con- 
vention. Our results are about 3 times larger than the 
smallest estimate amon g va rious QCD model predictions, 
|a| = \f3\ = 0.003 GeV3[l3|, which is often used in phe- 
nomenology of GUTs to derive "conservative" estimates 
of proton lifetime. 

Our present a and /? are smaller than the previ- 
ous results using the same gauge and quark actions, 
|a(NDR,l/a)| = 0.015(1) GeV^ and l^(NDR,l/a)| = 
0.014(1) GeV3 at 1/a = 2.30(4) GeV 0, beyond what 
is expected from scaling violation obtained in this work. 
We suspect that ^/Z^ and the lattice scale determined 
from nip are overestimated while their errors are not 
properly estimated, probably due to large fluctuations 
in ^yz^ for which only 100 configurations were used. 
We also compare our a and /3 with those of the pre- 
liminary results of RBC collaboration using quenched 
domain wall QCD with the DBW2 gauge action on an 
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FIG. 3: Continuum extrapolation of a (open) and /? (filled) 
parameters. The lattice scale is determined from rriM (circle), 
nip (square) and /,r (triangle) at each f3. The errors in the 
continuum limit are statistical only. 



L^xTxTVs = 16^ X 32 X 12 lattice at l/a= 1.23(5) GeV: 
|a(NDR,l/a)| = 0.006(1) GeV^ and |/3(NDR, l/a)| = 
0.007(1) GeV^ li, which are smaller by 30% than our 
values. This is of the order of scaling violation that is 
expected when the physical scale is set by mesonic quan- 
tities, but a quantitative comparison awaits their calcu- 
lation of the continuum limit. It should be noted that 
the effect of the change of renormalization scale from 
1/a = 1.23 GeV to 2 GeV is about 3.5% and negligi- 
bly small. 

In conclusion, we have studied scaling violation in the 
proton decay a and (3 parameters in quenched lattice 
QCD. Scaling violation is mild if the physical lattice scale 
is set by the nucleon mass, whereas a 30% systematic er- 
rors may arise from the physical scale, reflecting a part of 
the quenching error. Our estimate of a and /? is larger by 
about 3 times than the smallest prediction among various 
QCD models. This implies stronger constraints on GUT 
models. We can eliminate a remaining major uncertainty 
due to the quenched approximation by repeating the cal- 
culation on the full QCD gauge configurations which we 
already have|l^. This should be a next task. 
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